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The  operating  temperature  of  sodium-sulfur  battery  cells  is  above  300  °C  for  use  in  molten  liquid  state 
electrodes.  In  order  to  achieve  this  high  operating  temperature  condition  in  a  battery  module  composed 
of  multiple  cells,  a  thermal  management  system,  such  as  an  electric  heater  and  insulation,  is  essential. 
The  efficiency  of  the  module  is  directly  dependent  upon  the  temperature  uniformity  inside  the  module 
and  the  heat  dissipation  from  the  casing.  In  the  present  study,  a  new  numerical  model  for  the  thermal 
analysis  of  a  sodium-sulfur  battery  module  is  suggested.  The  equivalent  thermal  properties  of  the  cell 
are  evaluated  by  detailed  thermal  analysis  on  the  cell.  The  heat  generation  of  the  cell  is  modeled  con¬ 
sidering  the  electrochemical  reaction  process  and  the  variation  of  the  resistance  of  the  battery.  Using 
these  equivalent  thermal  models  of  the  cell,  a  zero  dimensional  lumped  thermal  model  for  examining 
the  effects  of  insulation  and  heater  operation  is  developed.  Finally,  the  three-dimensional  temperature 
distribution  inside  the  battery  module  is  predicted  by  solving  the  thermal  energy  conservation  equation 
numerically.  The  distribution  of  temperature  and  the  thermal  energy  efficiency  of  the  battery  module  for 
various  design  variables,  such  as  cell  arrangement  and  heater  operations,  are  summarized. 

©  2012  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

The  sodium-sulfur  battery  is  widely  known  for  having  a  high 
energy  density,  high  charge/discharge  efficiency,  and  long  cycle 
life.  Since  the  fundamental  research  on  this  battery  was  carried  out 
by  Ford  Motors  [1]  in  the  1960s,  the  early  studies  on  this  battery 
were  focused  on  exploiting  these  properties  for  application  to  elec¬ 
tric  vehicles  in  Europe  and  United  States  [2].  Subsequently,  since 
the  1980s,  the  battery  has  also  been  considered  a  promising  candi¬ 
date  for  stationary  energy  storage.  Indeed,  great  achievements  have 
been  made  for  this  application,  especially  in  Japan  [3].  Recently, 
researchers  in  China  [4]  are  also  carrying  out  the  study  on  the 
beta-alumina  ceramic  tubes. 

The  basic  working  principle  of  this  battery  is  the  electrochemical 
reaction  between  the  molten  sodium  (cathode)  and  sulfur  (anode) 
electrodes  [5]: 

discharge 

2Na  +  xS  pi  Na2Sx  (1) 

charge 

To  keep  sodium  and  sulfur  in  the  liquid  state,  the  cell  must  be 
operated  at  high  temperatures,  that  is,  in  the  range  290-350  °C.  At 
the  higher  end  of  this  range  (350  °C),  the  electromotive  force  (EMF) 
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of  the  cell  reaches  2.076  V  for  the  two-phase  region  and  declines  to 
1.78  V  as  the  depth  of  discharge  increases. 

For  a  practical  device,  to  achieve  the  necessary  operation  con¬ 
ditions,  each  battery  “module,”  comprising  multiple  unitary  cells, 
must  be  provided  with  a  thermal  management  system  including 
devices  such  as  heaters  and  an  insulating  casing.  According  to  Eck 
[6],  this  thermal  management  system  should  be  able  to  heat  the  bat¬ 
tery  to  the  desired  temperature  and  maintain  an  even  temperature 
distribution  under  all  operating  conditions— that  is,  the  maximum 
temperature  difference  should  not  exceed  25  and  20  °C  in  the  hor¬ 
izontal  and  vertical  directions,  respectively.  The  required  size  of 
the  battery  module  for  current  energy  storage  applications,  how¬ 
ever,  continues  to  increase,  making  it  increasingly  more  difficult 
to  maintain  temperature  uniformity  inside  the  module.  Another 
important  factor  is  the  heat  dissipation  from  the  module,  which  is 
directly  related  to  the  module  efficiency  [7].  As  a  result,  optimal 
design  of  the  thermal  management  system  for  the  battery  module 
is  essential  for  practical  applications  of  this  battery  system. 

To  this  end,  numerical  models  that  can  predict  the  thermal  per¬ 
formance  of  a  battery  module  would  be  useful  tools  in  the  design 
of  a  thermal  management  system.  Gu  and  Wang  [8]  summarized 
the  general  form  of  a  thermal-electrochemical  model  for  a  bat¬ 
tery  system  and  applied  it  to  a  Ni-MH  battery.  Such  models  have 
also  been  widely  used  for  fuel  cell  applications  [9-11].  Indeed, 
several  commercial  computational  fluid  dynamics  (CFD)  software 
packages  now  include  a  specialized  library  module  for  fuel  cell 
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Table  1 

Specifications  of  the  unitary  sodium-sulfur  cell. 


Specifications 

Note 

Capacity  [Ah] 

740.0 

Dimension  [mm] 

D  =  90.5/H=  542 

Resistance  [m£2] 

1.1 

Current  [A] 

(D/C)  80.6/(C)  71.4 

Operation 

88-650  Ah 

DoD  12-88% 

thermal-electrochemical  analysis  [12].  For  the  sodium-sulfur  bat¬ 
tery,  in  particular,  Hussein  et  al.  [13]  carried  out  a  simulation  to 
predict  the  voltage-current  behavior  and  also  included  a  simple 
model  of  the  temperature  variation  of  the  module  by  considering 
the  internal  resistance  variation. 

Although  there  has  been  great  progress  in  computer  hardware 
performance,  the  current  technology  is  not  capable  of  including 
the  detailed  shape  of  each  unitary  cell  and  its  corresponding  heat 
transfer  characteristics  in  the  module  calculation.  The  conventional 
alternative  is  to  use  the  equivalent  thermal  properties  [14].  This 
is  a  useful  concept,  especially  for  thermal  analysis  of  large-scale 
geometries  that  include  many  periodic  internal  substructures  with 
a  complex  configuration.  Another  important  aspect  of  engineering 
design  with  numerical  simulation  is  that  the  computational  pro¬ 
cedure  proceeds  in  a  stepwise  fashion  with  several  multi-fidelity 
sub-models.  At  the  early  stages  of  the  design,  the  designer  can 
narrow  down  the  important  design  variables  by  using  low-fidelity 
solvers  that  can  analyze  many  cases  in  a  relatively  short  period.  As 
a  result,  the  number  of  analyses  required  for  the  time-consuming 
detailed  high-fidelity  solver  can  be  minimized.  This  multi-step  and 
multi-fidelity  approach  is  the  key  strategy  for  systematic  engineer¬ 
ing  design. 

This  paper  presents  such  a  numerical  prediction  model  for 
sodium-sulfur  battery  modules.  The  equivalent  thermal  proper¬ 
ties  for  a  unitary  cell  are  evaluated  by  solving  the  steady-state 
heat  transfer  problem.  Then,  for  the  entire  module,  three  multi- 
step  computational  models  are  proposed.  The  low-fidelity  models 
adopt  lumped  capacitance  assumptions,  and  the  final  high-fidelity 
model  can  predict  the  detailed  thermal  performance  of  the  module. 
Thermal  analysis  on  a  typical  sodium-sulfur  battery  module  was 
carried  out  using  the  proposed  method.  The  effects  of  the  insulation 
wall  and  heater  operation,  module  efficiency,  and  the  temperature 
distribution  inside  the  module  are  summarized  quantitatively. 


2.  Development  of  numerical  models 

Fig.  1  represents  the  configurations  of  the  sodium-sulfur  cell 
and  module  considered  in  this  study.  Basically  the  model  was 
taken  from  the  open  literature  [15,16]  presented  by  NGK  Insula¬ 
tors.  Fig.  1(a)  shows  the  size  and  internal  shape  of  the  unitary  cell, 
which  contains  the  sodium  electrode  in  the  core,  the  electrolyte 
made  of  beta-alumina,  and  the  sulfur  electrode  surrounding  the 
electrolyte.  It  was  assumed  that  the  sodium  cartridge  is  made  of 
stainless  steel  (STS430)  and  the  material  of  the  cell  container  is  alu¬ 
minum  (AL3003).  By  connecting  these  cells,  the  module,  which  has 
a  capacity  of  360 kWh,  was  modeled  as  shown  in  Fig.  1(b),  which 
represents  the  in-line  cell  arrangement  at  the  horizontal  mid-plane 
of  the  module.  The  casing  of  the  module  is  made  of  insulation  mate¬ 
rial.  It  is  assumed  that  the  gaps  between  the  cells  are  completely 
filled  with  engineering  sand  made  of  SiCV  The  detailed  specifi¬ 
cations  of  the  cell  and  module  used  in  the  present  analysis  are 
summarized  in  Tables  1  and  2,  respectively.  For  the  information  not 
specified  in  Refs.  [15,16],  commercially  available  data  are  assumed. 
For  the  discharge  and  charge  period,  constant  current  operation  is 
assumed. 


Fig.  1.  Configurations  of  the  sodium-sulfur  battery  cell  and  module  considered  in 
this  study:  (a)  the  half-cell  configuration  and  (b)  the  horizontal  mid-plane  of  the 
module. 


As  a  systematic  numerical  approach  to  predicting  the  per¬ 
formance  of  a  thermal  management  system,  “multi-step  and 
multi-fidelity”  numerical  models  are  proposed  in  this  paper.  There 
are  a  cell  model  and  a  module  model  for  the  thermal  analysis. 
The  module  model  has  three  steps,  (i)  the  steady-state  lumped 
model,  (ii)  the  transient  lumped  model,  and  (iii)  the  detailed  three- 
dimensional  transient  model.  For  the  transient  thermal  analysis,  an 
electrochemical  model  for  the  evaluation  of  cell  heat  generation  is 
used. 

2.1.  Cell  model 

The  purpose  of  the  cell  analysis  is  to  evaluate  the  equivalent 
thermal  properties  of  the  unitary  cell:  the  thermal  conductivity, 


Table  2 

Specifications  of  the  sodium-sulfur  battery  module. 


Specifications 

Note 

Configuration 

8s  x  5p  x  8b  =  320  cells 

In-line 

arrangement 

Dimension  [mm] 

2200  x  1760x650 

Power  [kW] 

50 

Capacity  [kWh] 

360 

Operation 

7  h  (D/C)-5  h  (I)-8  h  (C)-4  h  (I)  =  24  h 
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(a)  (b) 


Fig.  2.  Cell  model  for  the  equivalent  thermal  property  calculation:  (a)  solid  modeling 
and  (b)  the  mesh  configuration. 


specific  heat,  and  the  density.  A  detailed  three-dimensional  heat 
conduction  analysis  was  adopted  for  this  end.  Fig.  2(a)  and  (b) 
shows  the  solid  model  of  the  cell  and  its  corresponding  mesh  con¬ 
figuration.  Using  the  symmetry  of  the  geometry,  only  half  of  the 
cell  was  modeled  and  it  was  assumed  that  the  cell  was  surrounded 
by  sand.  The  sodium  inside  the  core  part  of  the  cell  was  assumed 
to  have  a  50%  level.  After  a  detailed  mesh  dependency  test,  a  total 
of  956 1<  mesh  cells  were  generated  in  the  computational  domain. 
Using  the  commercial  Computational  Fluid  Dynamics  code  ANSYS- 
Fluent  6.3,  the  following  steady-state  heat  conduction  equation  was 
solved. 


d_ 

dxt 


=  0 


(2) 


where  k  is  the  thermal  conductivity  of  the  material,  T  the  tem¬ 
perature,  and  Xi  the  spatial  coordinate  in  the  i-direction.  For  the 
numerical  calculation,  the  finite  volume  approach  using  a  second- 
order  central  difference  scheme  was  used. 

As  it  can  be  seen  from  the  cylindrical  geometry  of  the  unitary 
cell,  the  thermal  conductivity  of  the  cell  should  be  different  in  the 
vertical  and  horizontal  directions.  In  other  words,  the  thermal  con¬ 
ductivity  should  have  an  anisotropic  or  an  orthotropic  property. 
In  the  present  calculation,  it  is  assumed  that  there  is  an  artificial 
temperature  difference  in  the  vertical  and  horizontal  directions, 
respectively.  The  equivalent  thermal  conductivity  of  the  computa¬ 
tional  model  in  the  i-direction  /<modei,i  can  be  evaluated  using  the 
following  Fourier’s  law  of  conduction. 

AT 

Qi  =  'WleuA^  (3) 

where  Qj  is  the  heat  transfer  rate,  the  cross  sectional  area,  A Tz 
the  temperature  difference,  and  Axz  the  length  of  the  model  in 
the  i-direction,  respectively.  The  conductivity  of  the  model  can  be 
expressed  using  the  volume  fraction  of  the  cell  or. 


^model  ,i  —  (1  °0^sand  +  °^cell,i 


(4) 


^  _  V°^ce\\ 

Vo  /model 


(5) 


Flere,  Volce\\  and  Vo/model  represent  the  volume  of  the  cell  and 
the  total  computational  domain,  respectively.  Using  Eq.  (4),  the 
equivalent  thermal  conductivity  of  the  cell  can  be  calculated  by 


Fig.  3.  Schematic  of  the  lumped  thermal  model  for  module  analysis. 


applying  the  finite  temperature  difference  in  the  vertical  and  hor¬ 
izontal  directions  as  the  boundary  condition.  At  the  other  sides  of 
the  domain,  the  thermally  adiabatic  condition  is  applied.  The  other 
equivalent  thermal  properties,  such  as  the  density  and  the  specific 
heat  of  the  cell,  can  be  calculated  by  the  volume  average  of  each 
part  of  the  cell. 

2.2.  Steady-state  lumped  thermal  model 

The  determination  of  insulation  material  and  its  geometric  spec¬ 
ifications  such  as  the  thickness  is  one  of  the  most  important  steps 
in  the  design  of  the  thermal  management  system  for  the  battery 
module.  Assuming  the  module  temperature  is  constant  inside  the 
casing,  as  shown  in  Fig.  3,  the  heat  dissipation  from  the  module  on 
the  z-th  side  can  be  calculated  by  solving  the  following  equation. 

^  _  \  ^  .  a  ^module  —  ^surface,! 

^-^dissipation  —  /  J  ^x- 

—  ^  ^  Mi(7surface,i  —  ^ambient)  (6) 

where  /q  is  the  thermal  conductivity  of  the  insulation  material, 
Ax*  the  thickness  of  the  insulation,  hi  the  convective  heat  transfer 
coefficient  of  z-th  side  of  the  module,  respectively.  A  simple  calcu¬ 
lation  program  has  been  generated  for  this  prediction.  By  changing 
the  thermal  conductivity  of  the  insulation,  the  thickness  of  insula¬ 
tion,  and  the  ambient  cooling  conditions,  the  total  amount  of  heat 
transfer  and  its  distribution  on  the  surfaces  of  the  module  can  be 
predicted. 

2.3.  Electrochemistry  of  the  battery  cell 

For  the  transient  thermal  analysis  of  the  module,  the  heat  gener¬ 
ation  inside  each  cell  should  be  taken  into  consideration.  The  heat 
generation  inside  the  cell  can  be  represented  by  following  equation 
[17-19]. 

Qcell  —  Qjoule  +  Qreaction  —  I  (^1  ~  (2) 

where  I  is  the  current,  r\  the  polarization  and  ( dE/dT )  is  the  temper¬ 
ature  dependence  of  the  open  circuit  voltage.  The  first  term  of  Eq. 
(7)  represents  the  joule  heating  due  to  the  electric  resistance  of  the 
cell,  and  the  polarization  is  written  by 

rj  =  IR 


(8) 
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Depth  of  discharge  [%  w.r.t.  Na2S3] 


Fig.  4.  Variation  of  the  entropy  term  of  the  sodium-sulfur  cell  reaction. 


Qjoule  =  I2R  O) 

where  R  is  the  electric  resistance  of  the  cell.  During  the  discharge 
and  charge  process,  this  resistance  value  varies  according  to  the 
depth  of  discharge  (DoD).  For  large  capacity  battery,  this  is  mainly 
caused  by  the  variation  of  the  reaction  surface  inside  the  cell.  For 
poorly  designed  capillary  wick  tube  inside  the  cell,  the  cell  resis¬ 
tance  increases  as  the  area  of  the  reaction  surface  on  the  electrolyte 
decreases.  The  resistance  variation  can  be  obtained  from  the  uni¬ 
tary  cell  experiment  measuring  the  voltage  at  a  constant  current 
using  the  following  equation: 


R  = 


OCV-V 

1 


(10) 


where  OCV is  the  open  circuit  voltage  and  V is  the  measured  voltage. 
In  the  present  study,  the  variation  of  resistance,  which  is  calculated 
by  Eq.  (10),  is  used  in  the  calculation. 

The  second  term  in  Eq.  (7)  is  the  heat  of  reaction  term  and 
T(dE/dT )  is  called  the  entropy  term.  There  have  been  several  exper¬ 
imental  approaches  to  evaluate  this  entropy  term  using  adiabatic 
discharge/charge  measurement  [5,1 7,18].  Fig.  4  represents  the  typ¬ 
ical  variation  of  the  entropy  term  taken  from  previous  studies.  At 
DoD  =  60%,  there  is  a  sudden  change  of  the  pattern  and  this  is  due 
to  the  transition  between  the  two-phase  and  single-phase  regions. 
It  should  be  noted  that  this  sodium-sulfur  reaction  is  exothermic 
for  the  discharge  process  and  endothermic  for  the  charge  process, 
respectively.  In  the  present  calculation,  the  reaction  model  sug¬ 
gested  by  Knoedler  [17]  was  adopted  for  the  transient  thermal 
analysis. 


2.4.  Transient  lumped  thermal  model 

It  is  well  known  that  the  lumped  heat  transfer  model  always 
gives  efficient  results  for  the  thermal  design  [12]  as  shown  in  the 
steady-state  lumped  thermal  model.  By  extending  the  concept  to 
the  transient  case,  the  average  temperature  variation  of  the  mod¬ 
ule  according  to  the  discharge  and  charge  operation  schedule  can 
be  predicted  efficiently.  This  calculation  can  evaluate  the  module 
efficiency  considering  the  heat  dissipation  loss  as  well. 

Assuming  a  uniform  temperature  inside  the  module,  the  gov¬ 
erning  equation  for  the  transient  heat  balance  of  the  module  can 
be  written  as 

pCpVol—  =  Qheater  +  Qjoule  T  Qreaction  —  Qconvection 

(  dE\  (11) 

—  Qheater +1  (  ^  )  _  ^  A  A  (^surface- ^ambient) 


where  p  and  Cp  are  the  equivalent  density  and  specific  heat  of  the 
module,  respectively,  and  Qheater  is  the  amount  of  heat  by  elec¬ 
tric  heater.  To  solve  this  zero-dimensional  ordinary  differential 
equation,  a  numerical  code  was  developed  using  an  explicit  time 
integration  scheme.  The  present  tested  battery  module  has  24  h 
operation  schedule  which  is  shown  in  Table  2.  The  average  tem¬ 
perature  of  the  module  and  its  corresponding  heat  generation  and 
dissipation  were  predicted  over  this  period. 

2.5.  Three-dimensional  thermal  model 

The  final  step  of  the  calculation  is  the  three-dimensional 
transient  thermal  model  that  predicts  the  detailed  temperature 
distribution  inside  the  module  considering  a  three-dimensional 
geometric  effect. 

Assuming  the  gap  in  the  module  is  completely  filled  with  sand, 
the  governing  equation  for  the  conduction  analysis  is  written  as 

(,2) 

where  q is  the  volumetric  heat  source  term.  When  there  is  empty 
space  inside  the  module  that  is  not  filled  with  sand,  there  should 
be  natural  convection  of  air  and  the  transient  computation  with 
convection  becomes  extremely  time-consuming.  In  the  case,  how¬ 
ever,  the  concept  of  effective  thermal  conductivity  [20]  for  the 
convection  region  would  be  applicable  to  the  transient  conduction 
analysis.  The  same  numerical  scheme  used  for  the  cell  model  is  also 
used  in  the  analysis  and  the  implicit  time  integration  is  carried  out 
numerically.  To  impose  the  Joule  heating,  the  cell  reaction  and  the 
electric  heating  in  the  heat  source  term  for  each  control  volume  of 
the  generated  mesh,  a  numerical  code  was  developed  to  customize 
the  commercial  code  Fluent  6.3  using  the  User  Defined  Function 
(UDF)  capability  of  the  program. 

Fig.  5(a)  shows  the  solid  model  of  the  current  tested  module. 
Using  the  symmetry  of  the  geometry,  only  a  quarter  part  of  the 
module  is  modeled.  Fig.  5(b)  represents  the  corresponding  mesh 
configuration  of  the  quarter  model,  which  has  1.2  M  hexahedral 
control  volume.  In  the  model,  the  cells  are  assumed  to  have  the 
equivalent  thermal  properties  evaluated  in  the  cell  model  calcula¬ 
tion. 

In  all  the  calculations,  the  ambient  temperature  is  set  to 
^ambient =  20  °C  and  the  convection  heat  transfer  coefficient  is 
h  =  5Wm-2K-1.  The  internal  temperature  of  the  steady-state 
lumped  model  Tmoduie  was  assumed  to  have  350  °C,  and  the  ini¬ 
tial  temperature  of  the  transient  lumped  model  Tmodule  was  set  to 
310  °C.  Since  the  current  thermal  models  consider  the  conduction 
heat  transfer  only,  and  the  internal  heat  generation  of  the  cell  due 
to  the  electric  resistance  and  electrochemical  reaction  is  very  slow, 
the  time  step  for  the  present  numerical  integration  could  be  rela¬ 
tively  large  compared  to  the  heat  transfer  problem  with  convection. 
Szabo  [21  ]  summarized  the  necessary  time  step  size  for  numerical 
stability.  In  the  present  transient  calculations,  a  time  step  of  A  t  =  5  s 
was  enough  to  satisfy  the  stability  criterion  and  independency  of 
the  time  step  size. 

3.  Results  and  discussion 

3. 1 .  Equivalent  thermal  property  of  the  cell 

Fig.  6(a)  and  (b)  shows  the  calculated  temperature  distributions 
inside  the  cell  for  vertical  and  horizontal  temperature  loadings, 
respectively.  The  temperature  difference  in  each  direction  was  set 
to  100°C.  Using  the  method  described  in  the  previous  section  for 
the  cell  model,  the  orthotropic  thermal  conductivity  was  eval¬ 
uated.  Fig.  7  represents  the  calculated  thermal  conductivities  in 
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Fig.  5.  Module  model  for  the  detailed  3D  calculation:  (a)  solid  modeling  and  (b) 
mesh  configuration. 


the  vertical  and  horizontal  directions  with  different  casing  mate¬ 
rials.  Here,  the  conductivity  is  normalized  by  the  volume  average 
value  of  the  thermal  conductivity  of  the  cell.  The  normalized  val¬ 
ues  are  much  smaller  than  unity,  which  represents  the  importance 
of  the  current  cell  model  approach  to  get  the  direction  dependent 
thermal  conductivity.  The  conductivity  in  the  vertical  direction 
shows  a  larger  value  than  that  of  the  horizontal  direction.  For  a 
stainless  steel  casing  it  is  larger  by  170%  and  for  aluminum  cas¬ 
ing  by  123%.  This  is  due  to  the  cylindrical  shape  of  the  internal 
structure.  The  shell  type  internal  casing  and  electrolyte,  which  are 
aligned  parallel  in  the  vertical  direction,  deliver  the  heat  along  the 
longitudinal  direction.  In  the  horizontal  direction,  however,  each 
component  is  separated  by  the  others,  which  results  in  a  larger 
thermal  resistance  connected  in  series.  When  the  casing  material 
is  changed  to  aluminum,  there  is  an  increase  of  conductivity  by  36% 
and  64%  in  the  horizontal  and  vertical  directions,  respectively.  The 
other  thermal  properties,  such  as  density  and  specific  heat,  can  be 


Fig.  6.  Temperature  distribution  inside  the  cell  under:  (a)  vertical  temperature  load¬ 
ing  and  (b)  horizontal  temperature  loading. 


calculated  by  taking  the  volume  average  values  and  using  the  vol¬ 
ume  fraction  relation.  The  baseline  cell  model  made  of  aluminum 
casing  in  the  present  study  has  an  equivalent  thermal  conductivity 
of  5.39  W  nrr 1  K-1  and  2.42  W  nrr 1  K-1  in  the  vertical  and  horizon¬ 
tal  directions,  respectively.  The  equivalent  density  and  the  specific 
heat  are  2247.1  kg  nrr3  and  1103.8Jkg-1  K_1,  respectively. 

3.2.  Heat  balance  and  insulation 

NGK  insulators  [15]  reported  that  the  insulation  of  the  top  face  of 
the  module  should  be  thinner  than  other  faces  in  order  to  reduce  the 
temperature  difference  in  the  horizontal  direction.  In  the  current 
calculation  using  a  steady-state  lumped  model,  two  different  insu¬ 
lation  thicknesses  were  tested:  (i)  50  mm  for  everywhere  except 
the  top  face  which  has  4  mm  (baseline  model)  and  (ii)  1 5  mm  every¬ 
where.  For  the  insulation  material,  the  typical  thermal  conductivity 
of  vacuum  insulation  was  used  (/<  =  0.01  Wm_1  K_1).  These  two 
types  of  insulation  thickness  were  selected  based  on  the  heat  dis¬ 
sipation  capacity  described  in  NGK  reference  [15],  which  has  the 
range  of  2200-3200  W.  In  the  present  paper,  the  insulation  thick¬ 
nesses  were  selected  to  have  the  dissipation  of  2700  W  (the  mean 
value)  by  several  iterations.  Table  3  indicates  the  distribution  of 
heat  dissipation  from  the  module.  The  heat  dissipation  through 
each  face  is  directly  proportional  to  the  thermal  conductivity  of 
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Fig.  7.  Dimensionless  equivalent  thermal  conductivity  of  the  cell. 
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Table  3 

Distribution  of  heat  dissipation  from  the  module  for  different  insulation  thicknesses 
[%]. 


Faces 

Top  =  4  mm 

Else  =  50  mm 

1 5  mm  everywhere 

Top 

78.8 

30.0 

Bottom 

9.1 

30.0 

Left 

2.7 

8.9 

Right 

2.7 

8.9 

Front 

3.4 

11.1 

Back 

3.4 

11.1 

the  insulation,  area  of  the  face,  and  the  reciprocal  of  insulation 
thickness  at  a  given  external  cooling  and  internal  temperature.  As 
represented  in  the  table,  for  a  thin  top  wall,  78.8%  of  the  heat  is 
dissipated  from  the  top  surface  as  it  was  intended.  If  the  module 
used  uniform  insulation  thickness  for  all  the  face,  the  distribution 
directly  follows  the  area  ratio  as  is  shown  in  the  second  column  of 
the  table,  while  the  total  amount  of  heat  transfer  rate  of  2700  W 
is  almost  the  same  as  that  for  a  thin  top  wall.  The  effect  of  the 
insulation  material  property  on  the  total  heat  dissipation  was  also 
examined.  If  the  thermal  conductivity  of  the  insulation  wall  was 
changed  to  0.023  Wm-1  K_1,  the  total  heat  dissipation  became 
4670  W.  This  prediction  is  helpful  for  designing  heater  capacity  at 
the  early  stage  of  the  module  design. 

3.3.  Temperature  variation  and  heater  operation 

In  the  transient  lumped  thermal  model  for  the  battery  mod¬ 
ule,  the  ideal  heater  operation  can  be  evaluated  by  considering 
the  heat  balance  between  heat  generation  in  the  cell  and  the  dis¬ 
sipation  from  the  module.  Table  4  indicates  this  ideal  minimum 
heater  operation  of  the  tested  module.  Fig.  8(a)  shows  the  amount 
of  heat  generation  and  dissipation  of  the  baseline  module,  and 
the  variation  of  module  temperature  for  a  24  h  cycle  is  shown  in 
Fig.  8(b).  During  the  discharge  period,  the  amount  of  heat  gen¬ 
eration  due  to  the  Joule  heating  and  reaction  heating  inside  the 
cells  is  larger  than  that  of  heat  dissipation  from  the  casing.  This 
excessive  internal  heating  results  in  an  increase  of  module  tem¬ 
perature  and  there  is  no  need  for  additional  heating  by  electric 
heaters.  In  the  first  idle  period  (Idle#l),  the  elevated  tempera¬ 
ture  decreases  by  the  external  cooling.  At  the  end  of  the  first 
idle  process,  however,  the  average  temperature  of  the  module  is 
still  higher  than  the  initial  temperature,  so  the  heater  should  be 
off  for  this  period  as  well.  In  the  charge  process,  although  the 
endothermic  reaction  occurs  in  the  cell,  due  to  the  Joule  heating 
and  higher  initial  temperature,  the  temperature  does  not  decrease 
lower  than  the  initial  temperature.  Finally,  in  the  second  idle 
period  that  completes  the  24-h  cycle,  the  heater  operation  is 
necessary  to  prevent  a  temperature  decrease  below  the  initial  tem¬ 
perature.  The  exact  capacity  of  the  heater  can  be  calculated  using 
the  heat  balance  after  several  iterations. 

In  Fig.  8(b),  the  circular  symbol  represents  the  experimental  data 
from  Refs.  [15,16].  The  time  averaged  difference  of  temperature 
between  the  calculation  and  measurement  is  5  °C.  Considering  the 
limited  information  in  Ref.  [15],  the  agreement  between  the  two 
is  good.  The  detailed  information,  however,  such  as  the  environ¬ 
mental  conditions  of  the  module  under  operation,  the  location  of 
temperature  measurement  sensors,  the  capacity  of  the  heater  and 
its  operation  schedule,  and  the  thermal  properties  of  the  insulation 
wall  could  improve  the  agreement  better. 

Table  4  also  contains  the  heater  operation  when  the 
insulation  material  has  larger  thermal  conductivity  values 
(/<  =  0.023  Wm-1  K-1)  and  the  resulting  temperature  variation  is 
shown  in  Fig.  8(b)  with  the  dotted  line.  Due  to  relatively  poorer 
insulation  capability,  the  heat  dissipation  becomes  larger  than  in 


Fig.  8.  Results  of  the  transient  lumped  model  for  the  baseline  module 
(/<  =  0.01  W  m-1 1<-1 ):  (a)  the  heat  generation  and  dissipation  of  the  module  and  (b) 
variation  of  the  module  temperature  for  the  24  h  cycle. 


the  case  of  baseline.  As  a  result,  the  heater  operation  is  neces¬ 
sary  from  the  first  idle  period.  During  the  charge  process,  it  can 
be  seen  from  the  figure  that  there  is  a  region  in  which  the  temper¬ 
ature  drops  below  the  initial  temperature  due  to  the  endothermic 
reaction.  After  the  whole  cycle,  however,  the  temperature  recov¬ 
ers  the  initial  level  since  the  heater  capacity  was  set  to  balance  the 
dissipation  and  generation. 

3.4.  Efficiency  consideration  and  variable  insulation  thickness 

As  summarized  by  Okuyama  and  Nomura  [7],  the  efficiencies  of 
the  battery  and  the  module  are  represented  as  follows: 

»?i  =  g  (13) 

r]2  =  -fff.  (14) 

u  B+C  K  J 

where  A  is  the  discharge  energy  output,  B  the  charge  energy  input, 
and  C  the  heater  power  consumption.  The  efficiencies  of  the  tested 
models  in  this  study  are  also  represented  in  Table  4.  When  the 
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Fig.  9.  Temperature  distribution  inside  the  baseline  module  for  various  times  at  the  symmetric  planes,  the  vertical  mid-plane,  and  the  horizontal  mid-plane. 


insulation  has  lower  conductivity  (/<  =  0.01  W  m-1 1<-1 ),  the  module 
efficiency  is  close  to  the  battery  efficiency.  With  poorer  insulation 
(/<  =  0.023  W  m-1 1<-1 ),  however,  the  module  efficiency  decreases  by 
9%,  which  was  caused  by  the  additional  heater  operation.  It  should 
be  noted  that  the  current  module  efficiency  is  the  ideal  (theoretical) 
efficiency  since  the  minimum  heater  operation  was  assumed  at  the 
given  condition.  In  real  operation  conditions,  the  module  efficiency 
should  usually  be  smaller  than  this  ideal  one  because  of  the  varia¬ 
tion  of  external  conditions.  The  active  control  of  heater  operation 
by  monitoring  the  module  temperature  is  essential  for  practical 
and  reliable  operation.  Nonetheless,  still  the  current  ideal  module 
efficiency  is  important  for  the  thermal  design  as  it  represents  the 
fundamental  performance  of  the  module  system. 

The  concept  of  using  variable  conductance  insulation  for  the  bat¬ 
tery  module  is  to  optimize  the  dissipation  of  heat  from  the  casing  by 
controlling  the  conductivity  of  the  material.  The  necessary  optimal 
heat  dissipations  from  the  casing  for  the  discharge  period  should 
be  larger  than  those  for  the  charge  and  idle  processes  because  of 
exothermic  internal  cell  heat  generation.  According  to  the  previous 
studies,  the  control  of  heat  dissipation  was  possible  by  changing  the 
pressure  level  inside  the  vacuum  insulation  panel  [22].  The  same 
effect  is  also  possible  by  changing  the  thickness  of  the  casing  which 
has  fixed  thermal  conductivity,  as  the  heat  transfer  rate  through  the 


insulation  wall  is  proportional  to  the  thermal  conductivity  and  the 
reciprocal  of  thickness.  In  the  present  study,  the  variable  thickness 
model  was  tested,  and  the  results  are  indicated  in  the  last  column 
of  Table  4.  In  the  calculation,  the  upper  thickness  was  set  to  4  mm 
for  the  discharge  process  and  set  to  50  mm  for  other  periods,  which 
has  the  same  meaning  as  the  conductivity  decreased  to  8%  of  the 
original  value  by  changing  the  vacuum  level.  The  result  shows  there 
is  a  great  improvement  in  the  module  efficiency,  i.e.,  the  efficiency 
of  78%  for  fixed  thickness  increases  to  85.8%.  The  corresponding 
temperature  variation  is  also  shown  in  Fig.  8(b)  with  a  dash-dot 
line.  It  is  obvious  that  the  case  of  the  variable  insulation  thickness 
shows  the  smallest  temperature  variation  of  the  module. 

3.5.  Temperature  uniformity  inside  the  module 

Fig.  9  shows  the  temperature  distribution  inside  the  module 
at  the  end  of  discharge  (t=7h),  the  end  of  the  first  idle  period 
( t  =  1 2  h),  and  the  end  of  the  charge  period  ( t  =  20  h).  The  tempera¬ 
ture  becomes  highest  at  t  =  7  h  due  to  the  excessive  internal  heating, 
and  decreases  as  the  cycle  goes  since  the  heat  dissipation  by  exter¬ 
nal  cooling  is  larger  than  the  internal  heat  generation  for  the  rest 
of  the  process.  The  thin  top  wall  allows  the  largest  heat  dissipa¬ 
tion  through  the  wall.  As  a  result,  the  temperature  gradient  in  the 


Table  4 

Module  and  ideal  heater  operations  for  different  insulations. 


Mode 

Time  [h] 

Heater  operation  [W] 

fc  =  0.01  Wm-1  K-1 
(baseline) 

k  =  0.023  Wm-1  K-1 
(fixed  thickness) 

k  =  0.023  Wm-1 1<-1 
(variable  thickness) 

Discharge 

7 

0.0 

0.0 

0.0 

Idle_l 

5 

0.0 

680.0 

0.0 

Charge 

8 

0.0 

3360.0 

442.0 

Idle_2 

4 

2000 

3660.0 

1426.0 

m  =A/B 

87.8% 

87.8% 

87.8% 

r)2  =AI(B  +  C) 

86.0% 

78.7% 

85.8% 
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Fig.  10.  Variation  of  the  average,  maximum,  and  minimum  temperatures  for  a  24  h 
module  cycle. 

vertical  direction  is  larger  than  that  in  the  horizontal  direction 
as  it  was  intended.  Since  the  equivalent  thermal  conductivity 
of  the  cell  is  larger  in  the  vertical  direction,  generating  larger 
dissipation  in  the  vertical  direction  is  a  better  strategy  for  mini¬ 
mizing  the  temperature  gradient  as  much  as  possible.  Due  to  the 
uneven  insulation  thickness,  the  maximum  temperature  occurs 
relatively  near  the  bottom  wall  in  the  central  core  region.  From 
the  result  of  t=7h,  it  can  be  clearly  seen  that  the  temperature 
of  the  cell  is  higher  than  the  surrounding  sand  region  due  to 
the  internal  heat  generation  by  resistance  and  battery  reaction. 
In  the  horizontal  direction,  the  internal  temperature  decreases 
from  the  center  to  the  side  wall  gradually.  For  the  cells  near 
the  corner  edge  region  of  the  module,  the  temperature  becomes 
lowest,  since  the  external  cooling  surface  area  for  those  cells  is 
larger  than  that  for  other  side  cells.  Consequently,  the  corner 
cells  have  the  worst  operating  conditions,  i.e.,  relatively  lower 
temperature  and  larger  temperature  gradient  in  both  directions. 
Since  the  temperature  of  the  corner  cell  is  largely  dependent 
on  the  distance  between  the  cell  and  the  insulation,  the  gap 
size  should  be  an  important  design  variable  and  carefully  deter¬ 
mined. 

Fig.  10  represents  the  variation  of  various  representative  tem¬ 
peratures  inside  the  module  for  a  24  h  cycle.  The  variation  of  the 


average  module  temperature  represented  by  a  solid  line  is  com¬ 
pared  with  the  prediction  results  by  the  transient  lumped  model 
shown  by  a  circular  symbol.  The  maximum  difference  between 
those  two  calculations  was  7.5  °C  and  the  average  difference  was 
3.0  °C.  The  average  cell  temperature  variation  is  close  to  that  of  the 
module;  it  is  slightly  higher  than  the  module  temperature  for  all 
cycle  as  the  cell  has  the  heat  generation  and  heat  capacity  to  store 
it.  Since  it  was  assumed  that  the  initial  temperature  was  uniform, 
the  difference  between  the  maximum  and  minimum  temperatures 
becomes  larger  as  the  cycle  proceeds.  Even  though  this  difference 
reaches  53  °C  at  the  end  of  the  charge  process  ( t  =  20  h),  it  should  be 
noted  that  these  maximum  and  minimum  values  were  measured  at 
the  local  extreme  points  inside  the  module.  For  example,  the  tem¬ 
perature  difference  of  the  cells  on  the  horizontal  mid-plane  is  less 
than  20  °C,  and  the  vertical  temperature  difference  on  the  vertical 
mid-plane  is  less  than  30  °C. 

In  the  present  result  for  the  baseline  model,  the  heating  by 
electric  heater  begins  from  t=20h  when  the  charge  process  ends. 
This  heating  prevents  the  module  temperature  from  decreasing 
lower  than  the  desired  level.  As  a  result,  the  module  and  the  cell 
temperature  do  not  drop  during  the  idle#2  period.  The  minimum 
temperature,  however,  increases  rapidly  since  the  corner  cell  is 
closely  located  near  the  heater,  which  represents  the  importance  of 
the  gap  size  near  the  side  and  corner  region  in  terms  of  temperature 
stability. 

Fig.  11(a)  shows  the  temperature  distribution  on  the  hori¬ 
zontal  mid-plane  when  the  conductivity  of  the  insulation  wall 
was  changed  to  a  larger  value  of  k  =  0.023  Wm-1  K_1.  As  shown 
in  Table  4,  for  this  case,  the  heater  heating  should  start  from 
the  charge  process  (t=12h)  due  to  the  large  heat  conduction 
through  the  casing  wall.  The  resulting  temperature  distribution 
at  t=20h  is  quite  different  from  the  result  of  t=12h.  By  electric 
heating,  the  temperature  in  the  side  and  corner  cells  increases, 
and  the  variation  of  temperature  on  the  horizontal  mid-plane 
becomes  smaller.  In  other  words,  the  use  of  electric  heating 
has  a  potential  to  improve  the  temperature  uniformity  with  the 
penalty  of  module  efficiency.  Fig.  11(b)  represents  the  result  of 
temperature  variation  on  the  horizontal  mid-plane  of  a  mod¬ 
ule  having  a  staggered  arrangement  of  cells.  When  the  cells 
are  arranged  in  staggered  way,  the  average  gap  between  the 
cells  decreases  and  more  compact  module  can  be  constructed. 
The  cooling  rate  of  the  module  is  directly  dependent  upon  the 
surface-to-volume  ratio.  When  the  volume  of  the  module  for  a 
fixed  number  of  cell  decreases,  the  volume  ratio  is  proportional 
to  the  cubic  power  of  the  length  scale,  whereas  the  module 
surface  decreases  proportional  to  the  square.  Consequently,  the 


t=7HR  (end  of  discharge) 
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Fig.  11.  Temperature  distribution  inside  the  module  for  different  cell  arrangements  with  an  insulation  of  k  =  0.023  Wirr1  K-1:  (a)  the  in-line  arrangement  and  (b)  the 
staggered  arrangement. 
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surface-to-volume  ratio  increases  linearly  to  the  length  ratio, 
which  means  the  cooling  speed  increases.  In  the  present  study 
shown  in  Fig.  11(b),  the  staggered  arrangement  has  smaller  vol¬ 
ume  by  7.9%  and  surface  area  by  6.6%.  As  a  result,  at  t=12h, 
more  cooling  effect  can  be  seen  in  the  figure  especially  near 
the  side  and  corner  regions.  In  terms  of  temperature  unifor¬ 
mity,  the  in-line  arrangement  has  an  advantage  when  the  heater 
is  operating,  as  is  shown  in  the  t  =  20h  results.  This  is  due  to 
the  relatively  slower  heat  conduction  speed  due  to  larger  ther¬ 
mal  resistance  and  heat  capacity  compared  to  the  staggered 
arrangement. 

4.  Conclusions 

A  fundamental  analysis  procedure  for  evaluating  the  ther¬ 
mal  performance  of  a  sodium-sulfur  battery  module  has  been 
developed  and  presented.  Due  to  the  complexity  of  the  mod¬ 
ule,  a  multi-step,  multi-fidelity  approach  is  essential  for  the 
systematic  design  process.  The  proposed  four-step  approach  is 
composed  of  (i)  3D  cell  analysis  for  evaluating  the  equiva¬ 
lent  thermal  properties  of  the  cell,  (ii)  a  steady-state  lumped 
thermal  model  for  estimating  the  heat  balance  and  insula¬ 
tion  thickness,  (iii)  a  transient  lumped  model  for  determining 
heater  operation  and  evaluating  module  efficiency,  and  (iv)  3D 
detailed  analysis  for  predicting  the  temperature  distribution  inside 
the  module.  For  the  cell  heat  generation  model,  the  varia¬ 
tion  of  the  resistance  during  the  discharge  and  charge  process 
was  taken  into  consideration,  and  the  electrochemical  reaction 
term  was  modeled  using  a  previous  correlation  for  the  entropy 
term. 

Using  the  method  suggested,  the  effects  of  various  design  vari¬ 
ables  were  examined  quantitatively.  The  material  property  of  the 
cell  and  the  corresponding  equivalent  property  of  the  cell  were 
evaluated  for  different  casing  materials.  The  orthotropic  equiva¬ 
lent  thermal  conductivity  of  the  cell  shows  larger  conductivity  in 
the  vertical  direction  than  in  the  horizontal  direction.  Distribution 
of  heat  dissipation  on  each  face  of  the  module  was  examined  for  dif¬ 
ferent  insulation  thicknesses  and  materials.  Module  efficiency  with 
various  heater  operations  and  insulation  materials  was  predicted. 
Using  the  variable  insulation  thickness  model,  the  improvement 
of  module  efficiency  is  clearly  predicted.  Temperature  distribu¬ 
tion  inside  the  module  for  different  operating  conditions  and  cell 
arrangements  was  predicted  using  detailed  3D  analysis.  The  impor¬ 
tance  of  gap  size  at  the  corner  and  side  regions  and  the  effect  of 
a  staggered  cell  arrangement  are  found.  The  variation  of  module 
temperature  for  initial  heating  period  to  elevate  the  tempera¬ 
ture  from  ambient  to  the  operating  temperature,  and  for  several 
periods  of  discharge-charge  cycles  will  be  incorporated  in  future 
work. 
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Appendix  A.  Nomenclature 

A  discharge  energy  output,  Wh 

Aj  cross  sectional  area  in  the  i-direction,  m2 

B  charge  energy  input,  Wh 

C  heater  power  consumption,  Wh 

Cp  specific  heat,  J  kg-1  K-1 

E  voltage,  V 

h  convective  heat  transfer  coefficient,  W  m-2  K-1 

/  current,  A 

k  thermal  conductivity,  W  m-1 1<-1 

Na  sodium 

OCV  open  circuit  voltage,  V 

qf//  volumetric  heat  source,  W  m-3 

Qi  heat  transfer  rate  in  the  i-direction,  W 

R  electric  resistance,  £2 

S  sulfur 

t  time,  s 

T  temperature,  K  or  °C 

V  measured  voltage,  V 

Vol  volume,  m3 

Xi  spatial  coordinate  in  the  i-direction,  m 

Greek 

a  volume  fraction 

At  time  step,  s 

AT/  temperature  difference  in  the  i-direction,  I<  or  °C 
Axj  thickness  of  the  model  in  the  i-direction,  m 
r)  polarization,  V 

r\\  battery  efficiency 

r\2  module  efficiency 

p  density,  kg  m-3 
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